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ABSTRACT 


The results of extensive computations are presented in order to accurately 
characterize transitions to chaos for the Kuramoto-Sivashinsky equation. In particular 
we follow the oscillatory dynamics in a window that supports a complete sequence of 
period doubling bifurcations preceding chaos. As many as thirteen period doublings 
are followed and used to compute the Feigenbaum number for the cascade and so 
enable, for the first time, an accurate numerical evaluation of the theory of universal 
behavior of nonlinear systems, for an infinite dimensional dynamical system. Further 
more, the dynamics at the threshold of chaos exhibit a fractal behavior which is 
demonstrated and used to compute a universal scaling factor that enables the self- 
similar continuation of the solution into a chaotic regime. 
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1. Introduction. 

A central question in fluid dynamics that is attracting a considerable research effort, is the predic- 
tion of onset to turbulence. A general theory encompassing the Navier-Stokes equations of fluid 
motion, and consequently covering a large class of physical phenomena, is not available at present. As 
a result most contributions are focused on the analysis of model equations derived from the Navier- 
Stokes system by asymptotic methods, for example, or by finite-dimensional truncations [1]. In many 
cases this is a valid and useful approach, especially in the light of Feigenbaum’s fascinating theory ori- 
ginally for one-dimensional nonlinear maps ([2]-[4]), which predicts universal nonlinear behavior and 
is believed to be applicable to many more complex nonlinear systems such as ordinary and partial 
differential equations. A brief review of Feigenbaum’s theory for the quadratic map is in order here, 
but the interested reader should refer to the above mentioned articles (also [5]). The theory pertains to 
a large class of mappings of an interval onto itself, a representative member of which has 

/(x) = 4vx(l-x) , v > 0 , x e [0,1] . (1) 

The flow map is represented by the repeated application of (1). This choice of (1) ensures that Xj = 0 
and x 2 = 1— l/4v are always fixed points (i.e. solutions of x = f (x)) and that for any 0 < V < 1 all 
iterates are contained in [0,1]. Now for 0 < v < 1/4 the fixed point x j alone is in the range of the 
map, while for 1/4 < v < 1 both fixed points Xj and x 2 are in the range. Stability of a fixed point is 
determined by the condition l/'(x)l < 1, which tells us that for 0 < V < 1/4 the fixed point Xj is 
stable, whereas for 1/4 < v < 3/4, x j is unstable and x 2 is stable (these results are easily obtainable 
by differentiation of (1)). Further more, at v = 3/4 we have f'(x 2 ) = -1 ( i.e. x 2 is a nort-hyperbolic 
point.) and so x 2 becomes unstable also for v > 3/4. To summarize, therefore, when 0 < V < 3/4 
repeated application of (1) starting from some initial position will yield convergence to one of the two 
fixed points and the flow is both deterministic and attractive to a global fixed point. Things become 
interesting when v is increased and lies in the range 3/4 < v < 1. For example, at values of v just 
above 3/4 a cycle which contains two alternating fixed points is obtained. In other words, repeated 
application of (1) gives a sequence x n , x 22 , x 12 , x 21 , • • • . Thus, a period doubling takes place 
and it is easy to find the values of Xj 2 and x 22 since they must be fixed points of the function 
f (f CO), also written as / 2 (x). Further increase of V makes these points unstable and a further split 
into a cycle containing 4 points heralds the second period doubling. Again the elements of this 4-cycle 
are computable by identification of the fixed points of / 4 (x). This process persists with period 
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doublings ad infinitum towards an accumulation point beyond which chaos sets in. More pre- 
cisely, if we define v„ to be the value of the parameter where the n* period doubling takes place, then 
the theory predicts that 

= v n+i ~ v " ... > lim 5„ = 5 = 4.6692016... . (2) 

V„ +2 - V n+ i » -» ~ 

The Feigenbaum number, 5, is universal for a variety of maps which have locally quadratic maxima. 
There is another important universal constant which arises from the theory of nonlinear maps. Without 
loss of generality, if we restrict our attention to the vicinity of the turning point x = 1/2 and follow a 
sequence of period doublings then it is found that the distance between neighboring elements of period 
doubled orbits (in the high-iteration limit) is reduced by a universal constant factor 

a = 2.502907875... . (3) 

For a comprehensive review of nonlinear maps see [6]. This property can be used as a predictive tool 
to construct the behavior of higher-iterates and even chaotic motions beyond the accumulation point of 
the cascade (see [7], for an application in Rayleigh-Benard flow). In what follows we provide con- 
clusive evidence that the universal behavior described above, also appears in the infinite-dimensional 
dynamical system provided by the Kuramoto-Sivashinsky equation, and both universal constants are 
computed with desirable accuracy. In achieving this, we obtain as many as 13 period doublings (the 
period increases from a value of for the first time as far as we know, conclusive evidence for a period 
doubling (classical) route to chaos for an infinite dimensional dynamical system. 

The equation studied is the Kuramoto-Sivashinsky equation (KS), suitably written in the form 

u, + uu x + u^ + V M yyTT = 0 , 

(x,t) e R ! xR + , (4) 

m(x,0) = m 0 (*) > u(x+2n,t) = u(x,t) , 

where v > 0 is the viscosity of the system. This equation arises in a variety of problems such as con r ? 
centration waves [8], flame propagation [9], free surface flows [10]. A generalized form, with (3 )^ls a 
special case, has been derived by an asymptotic analysis of the Navicr-Stokes in the context of two- 
phase flows in cylindrical geometries with applications in lubricated pipe-lining (for the efficient tran- 
sport of crude oil) and oil recovery through porous media [11]. Much analytical and computational 
work has been completed in order to elucidate some of the highly complex and nonlinear behavior that 
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(4) is capable of producing. For many references see [12] and [13]. The former article concentrated 
on a fairly global characterization of solutions as v varies, and in particular when it achieves fairly 
small values. Note that within the context of linear theory (for the periodic problem posed above), the 
number of unstable modes is equal to [V 1/2 ] ; as v decreases, therefore, more unstable modes enter 
and our concern is with the long-time behavior of these under the action of nonlinearity. 

2. Numerical solutions. 

The results presented here, were obtained by numerical solution of the initial value problem (4) 
with the initial condition 

h 0 (x) = -sin(x) , 

for all values of v. Inspection of (3) indicates that if the solution is an odd function of x initially, then 
it will remain so for all subsequent times. The advantage of such a choice is two-fold; first, the solu- 
tion can be expanded as a Fourier sine-series and a spectral Galerkin approximation is readily applica- 
ble - secondly, existing analytical results are available that give global bounds for u(x,t) (and higher 
derivatives) only in the odd-parity case [14], whereas analysis of the general case yields an estimate 
that the solution can grow at most exponentially [15]. The numerical scheme is described elsewhere as 
well as dynamics over a wide range of parameter values [12]. The truncation order of the Galerkin 
approximation depends on the value of v, but a crude estimate that has proven practical, is the reten- 
tion of a few frequencies more than the number of linearly unstable ones. This number is related to the 
dimension of the attractor for given v but no conclusive remarks can be made at this point (note that 
the theoretical upper bound estimate for the Hausdorff dimension of the attractor yields a value propor- 
tional to v- 21/40 which is larger than v~ 1/2 by a factor of v _1/40 , for details see [16] ; the multiplicative 
constant that comes with the bound, however, can be very large thus obstructing direct practical appli- 
cations). A brief summary of the phase-space is now given, and documented in Table 1. The descrip- 
tion we give starts at v=l and continues to ever-decreasing values. An initial bifurcation occurs when 
V decreases below 1 giving unimodal fixed points, and as v decreases further a series of Hopf bifurca- 
tions are encountered at v=l/4, 1/9,1/16,... which are seen as dimodal, trimodal, tetramodal etc. sta- 
tionary solutions. Not all orbits are stable, however, and so are not necessarily realized. What is more 
important is that instability appears that gives the solutions an oscillatory nature. A time-periodic 
attractor is obtained for the first time at a value of v ~ 0.06; a period-doubling occurs in this window 
also, but as v is decreased further the solutions are attracted to fully modal fixed points. A new time- 
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periodic window with two period-doublings and a period-halving (in this order), is then found 
v e [.03729, .03962] to separate the fully modal fixed points with another laminar window supporting 
tetramodal fixed points. The next periodic window v e [.029969 104,. 034625896] has the added 
feature that it contains a complete sequence of period-doubling bifurcations which lead to chaos. We 
call this the third periodic window and describe and analyze data from it presently. Note that the first 
chaotic window as well as the periodic one that preceded it, are only the first in a series as v is 
decreased. In fact at much smaller values of v we found more periodic windows with a complete 
sequence of period-doubling bifurcations leading to a chaotic flow. The lengths, in parameter space, of 
such successive attractors shrink significantly, but nonetheless we were able to compute the Feigen- 
baum number for three successive time-periodic attractors that lead to chaos. In what follows we 
describe the first periodic attractor that contains such a sequence, in order to make an accurate com- 
parison with Feigenbaum’s theory as applied to a partial differential equation. 

Com 2 utatwn of tjie Feigenbaum number. 

Table 2 presents our evidence that supports Feigenbaum’s universal theory for the Kuramoto- 
Sivashinsky equation. These results were generated by monitoring the time evolution of the energy, 
E(t), of the solution (the L 2 -norm) but any other positive definite quantity can be chosen also, for 
example the supremum of the solution. Each entry in Table 1 represents the beginning of successive 
sub-windows which support solutions that undergo period doublings. The sharp estimation of these 
boundaries is necessary if an accurate computation of the Feigenbaum number is desired. In all results 
provided here the boundaries were estimated with enough accuracy to yield the Feigenbaum number 
correct to the number of significant figures shown. The first column gives the value of v where the 
window begins, the second column gives the period of the oscillation and the third column gives the 
ratio of successive sub-window lengths according to equation (2). A more visual form of these results 
can be seen from a sequence of energy phase planes which cover the first 5 period doublings. It is 
clear from the figures that the overall bounds of the solution, for example the maximum and minimum 
of E and E , do not vary much beyond the second period doubling. Chaos sets in by the appearance of 
more turns in the phase plane, in other words by an index change of the curves (in fact the exact way 
in which the phase plane gains more turns before the appearance of chaos, is quantified in the next 
sub-section). Another qualitative aspect of the period doublings as seen through phase-plane 
sequences, is the index doubling of the curves at the beginning of successive sub-windows (the end of 



-5- 


super-stable limit cycles) by the splitting of one branch into two; the two branches then drift apart very 
slowly as time evolves until a new sub-window is encountered and another index doubling takes place. 
The time period itself increases very slightly over the duration of a sub-window but period doublings 
are exact and spontaneous. 

The universal limit of multiple period doublings^ . 

Next we present a set of numerical results that exhibit very clearly the fractal nature of the period 
doubling bifurcations. The experiment we choose has a value of v = .0299691035 and lies at the end 
of the third periodic window. The time period of the solution is 1798.2564595 units and is the result 
of a sequence of 12 period doublings (in Figure 1 we show only the first 5). Note that at a value of v 
=.029969103484, i.e. a decrease of 1.6xl0 -11 chaotic solutions were observed. Consider, then, the 
energy phase plane at v = .0299691035. This is not shown here but is similar to the last frame of 
Figure 1, except that it has more turns. We now look at that region of the phase plane which has 
!£■(/) I < 0.01 and whose horizontal axis is chosen so as to include initially all the minima of the 
energy time series (in other words all the branches on the left of the figure). This diagram is shown in 
Figure 2(a). Figure 2(b) is now obtained as follows : The minima on the right half of 2(a) are con- 
sidered and the picture is enlarged and centered in exactly the same way that 2(a) was (this ensures 
that the position of the leftmost and rightmost branches relative to a given enlarged picture are 
congruent). We now take the left half of 2(b) and enlarge it as indicated above to produce 2(c). The 
right half of 2(c) is then enlarged to produce 2(d) and so on. The reason for alternating between right 
and left halves is that we obtain the correct regions that need to be enlarged in order to demonstrate 
self-similarity of two symmetrically placed regions of the phase plane. This process provides a picture 
of the dynamics in a continuously shrinking region of the phase plane and is the appropriate procedure 
to follow in searching for any kind of self-similarity or fractal properties of the flow. The pictures 
presented here involve scrutinization of the "central" part of the phase plane (in fact keeping track of 
the sequence of magnifications we see that the portion of 2(a) which is magnified all the way to the 
end, is in the vicinity of E — 17). Other parts of the curve could have been chosen also but the 
present choice provided the simplest algorithm. With the above construction we are in a position to 
exhibit the fractal nature of the dynamics. The pictures 2(c), 2(e), 2(g) and 2(i) are seen to be self- 
similar as also are 2(b), 2(d), 2(f), 2(h) and 2(j). Self-similarity is deduced by the identical geometrical 
arrangement of the various portions of the phase plane as the magnification proceeds. Also, figures 
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which are side-by-side are reflections of each other about the vertical axis (this is to be expected by 
analogy with maps of quadratic nonlinearity). Next, consider the distance between the leftmost and 
rightmost elements of Figures 2(b)-2(j). If the ratio of successive distances is calculated as the 
magnification proceeds, the number found is 2.503 to within the accuracy of our measurements. This 
is the universal scaling factor a described in the Introduction (equation (3)), and so we have another 
instance of a complete confirmation of Feigenbaum’s universal theory for the Kuramoto-Sivashinsky 
equation. 
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Overview of the most attracting manifolds 

Window 

range 


Description of the attractors 

1 

< 

V 

< 

oo 

Constant states. 

.25 

< 

V 

< 

1 

Fully modal steady attractors. 

.0756 

< 

V 

< 

.25 

Fully modal and dimodal steady attractors. 

.06697 

< 

V 

< 

.0755 

Fully modal steady attractors. 

.05992 

< 

V 

< 

.06695 

Trimodal steady attractors. 

.05516 

< 

V 

< 

.05991 

Periodic attractors. 

.0396227 

< 

V 

< 

.05515 

Fully modal steady attractors. 

.03729 

< 

V 

< 

.0396226 

Periodic attractors. 

.0346259 

< 

V 

< 

.03728 

Tetramodal steady attractors. 

.029969103484 

< 

V 

< 

.0346258 

Periodic attractors. Complete period- 
doubling sequence. 

.02922 

< 

V 

< 

.02969910348 

Chaotic oscillations. 

.02905 

< 

V 

< 

.02921 

Periodic attractors. 

.02855 

< 

V 

< 

.02904 

Chaotic oscillations. 

.02662 

< 

V 

< 

.02854 

Periodic attractors. 

.02525 

< 

V 

< 

.02661 

Chaotic oscillations. 

.02506 

< 

V 

< 

.02524 

Periodic attractors. 

.0248607 

< 

V 

< 

.02505 

Chaotic oscillations. 

.02445 

< 

V 

< 

.0248606 

Periodic attractors. Complete period- 
doubling sequence. 

.0242861 

< 

V 

< 

.02445 

Chaotic oscillations. 

.02367 

< 

V 

< 

.02438608 

Periodic attractors. Complete period- 
doubling sequence. 

.0232 

< 

V 

< 

.02386 

Chaotic oscillations. 

.0229 

< 

V 

< 

.0231 

Periodic attractors. 

.0223 

< 

V 

< 

.0228 

Chaotic oscillations. 

.022 

< 

V 

< 

.0222 

Periodic attractors. 

? 

< 

V 

< 

.0219 

Chaotic oscillations. 


TABLE 1 
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20-mode Galerkin Expansion 

Subwindow boundary 

Length 

Ratio 

Period 

.0346258 

4.3083xl0~ 3 


.44 

.03031749 

2.6825X10 -4 


.88 

.030049233 

6.2786xl0~ 5 

16.061 

1.76 

.029986446 

1. 3609x1 r 5 

4.2724 

3.52 

.0299728366 

2.9330xl0" 6 

4.6136 

7.03 

.0299699036 

6.288xl0“ 7 

4.6399 

14.05 

.02996927484 

1.3456xl0“ 7 

4.6644 

28.1 

.02996914018 

2.884x1 0“ 8 

4.6657 

56.2 

.02996911134 

6. 18x1 0 -9 

4.667 

1 12.4 

.02996910516 

1.32xl0 -9 

4.68 

224.8 

.029969103842 

2.84xl0 _1 ° 

4.65 

449.6 

.029969103558 

6.0x1 0~ n 

4.7 

899.1 

.029969103498 

1.4xl0" n 

4. 

1798.2 

.029969103484 



oo 


TABLE 2 
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FIGURE 2 : 2(a)-2(j) - successive magnification of a portion of the phase-plane for the 2 12 -cycle, 
showing the fractal characteristics of the attractor. 
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